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Transcritical bifurcation of an 
immunosuppressive infection model 


E. Shamsara, R. Mostolizadeh and Z. Afsharnezhad* 


Abstract 


In this paper, the dynamic behavior of an immunosuppressive infection 
model, specifically AIDS, is analyzed. We show through a simple mathe- 
matical model that a sigmoidal CTL response can lead to the occurrence of 
transcritical bifurcation. This condition usually occurs in immunodeficiency 
virus infections (such as AIDS infection) in which viruses attack immune 
cells CD4*+T. Our results imply that the dynamic interactions between the 
CTL immune response and HIV infection are very complex and in the CTL 
response, dynamics can exist the stable regions and unstable regions. At the 
end of the paper, numerical simulations are presented to illustrate the main 
results. 
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1 Introduction 


One of the most complicated organs of higher organisms is the immune sys- 
tem. The function of the immune system is to fight off pathogenic organisms 
that enter and grow within the host (for example, viruses, bacteria, unicellu- 
lar eukaryotic parasites such as malaria, and multicellular parasites such as 
worms). Immune responses can be subdivided broadly into two categories: 
(i) innate or nonspecific responses, and (ii) specific, adaptive responses. In- 
nate immune mechanisms provide a first line of defense against an invading 
pathogen. They include physical barriers like the skin, changes in the envi- 
ronment of the body, such as fever, and immune cells that can fight pathogens 
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in a nonspecific way. Nonspecific is the key word here and means that these 
responses cannot specifically recognize the physical structure of the pathogen. 
Instead, these nonspecifics sense that an invader is present and react. While 
such responses slow down the initial growth of a pathogen, they are usually 
insufficient to clear an infection. For removing an infection, a specific and 
adaptive immune response is required [14]. The adaptive immune response 
consists of three main branches. 1. The B cells secret antibodies that neu- 
tralize free virus particles. 2. The CTL (also known as CD4* T cells) attack 
infected cells. 3. The CD4* T helper cells are very important regulators that 
ensure that CTL and B cell responses are developed efficiently. In immuno- 
suppressive infection models, infected cells attack to CD4* T cells and infect 
them; subsequently, they cannot help CTL and CD4* T cells to act efficiently. 
Mathematical models have been of central importance for understanding the 
dynamics between viral infections and immune responses, particularly in the 
context of a human immunodeficiency virus (HIV) infection [7]. Significant 
emphasis has been placed on the viral side of these dynamics, including the 
estimation of basic viral parameters. Subsequent work has focused on the 
immune side of these interactions in trying to explain a variety of experimen- 
tal observations about the dynamics of immune cells in various infections. 
One particular part of the immune system that is very important in the fight 
against viral infections is the killer T cells or cytotoxic T lymphocytes (CTL). 
They basically fight intracellular pathogens [15]. Clinical data have shown 
that for some human pathogens, such as HIV, hepatitis B virus (HBV) and 
hepatitis C virus (HCV), drug therapy sometimes is not completely effec- 
tive [7,13]. Recently, in 2015 [2,12] and 2014 [9], impaired immune responses 
in immunosuppressive infection models have attracted more and more atten- 
tion. Mathematical models have been developed to capture the interaction 
in vivo among HIV [2,4,5,9,12, 14, 15). 

The following model is general and satisfied the clinical data, so it was pur- 
sued by scientists; see the above references. In 2003, this model was developed 
and considered [5]. 


y = yor(y) — yz 
{ z= zf(y). @) 


In this system, y is the virus population and z is the population of the immune 
cells. The function g,(y) should be satisfied in: 


(2) 


1.g-(0) > 0, 92 <0Vy 
2.3y* >0,gr(y*) > 0, 2 > Ovr, y. 


Also the following conditions were assumed for f(y) in [5] 


3.4 y1, ye > Osuch that f(y) = f(y2) =0 
4. $f > 0 for y = yi and $f <0 for Y = Y2 
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gr(y) is the virus growth function that depends on the viral replication rate, 
r, and f(y) is the immune expansion function that does not depend on r. In 
the above case, when viral replication is high and the virus load is between 
y1 and yo, immune expansion is increasing and levels of antigen are sufficient 
to trigger sustained immunity [5]. In [5], a special function for g,(y) and f(y) 
is introduced (see Model (12)). Model (12), in 2015 [2,12], was considered to 
investigate the stability of the CTL immune response. Shu et al [9] in 2014 
obtained saddle point for system (12) which shows stable and unstable. Note 
that all the above investigations were on the eigenvalues with the non-zero 
real part and they didn’t consider the zero eigenvalue (bifurcation theory). 
In this paper, we are interested in only one zero eigenvalue of the system (12) 
at the fixed point, which can lead to the occurance of transcritical bifurcation 
[1,6, 8]. Since our concentration is on AIDS, we change the condition (3), 
in order to consider a weak immune system. The difference between HIV 
and AIDS is: HIV is the beginning of the AIDS disease, in AIDS; virus load 
rises more sharply, and the CD4* T cell (which defend against ADIS cells) 
drops sharply [14]. From a mathematical point of view, a > 0 means that 
the function f is a strictly increasing function with respect to the variable y. 
From a biological perspective, it means that the function of immune system 
responses to the disease increase. In this study the conditions for f(y) are as 
follows: 


L.Ay* >0;f(y*) =0 (4) 
2. $f = 0 for y* > 0. 


The new conditions cause a critical situation for the function f. For this case, 
we try to find a zero eigenvalue to apply transcritical bifurcation. Bifurcation 
theory helps us to obtain conditions for the parameters to keep the disease 
stable. In other words, by finding a region for parameter r with respect to 
parameter k, we tried to keep the immune system in proper condition as long 
as possible. Our work is organized as follows: 


In Section 2, we give some preliminary definitions of bifurcation and the- 
orems, which are going to be used in other Sections. Section 3 is devoted 
to bifurcation of system (12). Section 4 illustrates our numerical results. 
Section 5 is the conclusion. 


2 Preliminaries 


Bifurcation theory is fundamental for the qualitative study of dynamical sys- 
tems, and can be used to reveal complex dynamical behaviors of the biological 
systems under study, such as bistability, recurrence, and regular oscillation. 
Characterized by a controllable parameter, called the bifurcation parameter, 
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bifurcation occurs at a critical value of this parameter where the properties 
of equilibria change significantly. 


We consider bifurcations of equilibria of autonomous systems which de- 
pend on one single parameter p: 


&= f(z,p) ,2 ER” weR. (5) 


The system (5) is called smooth if f(x, 44) is differentiable up to any order in 
both # and yw. Equilibria of (5) are solutions of the algebraic equations 


f(x, w) = 0. (6) 


In order to graphically illustrate the dependence of an equilibrium x on pu, we 
require a scalar measure of the n-vector x. We shall use the notation [2] for 
such a measure of «. A diagram depicting [x] versus ju, where (2, ~) solves 
equation (1), will be called a bifurcation diagram. The continuous curves of 
solutions of (1) under variation of y are called branches. The branches of 
smooth systems are continuous and smooth but can split into more branches. 
On a regular point of a branch, that is, on a point where the branch does not 
split or turn around, we can define the slope of the branch. We will use the 
following abbreviations: 


Of (x, MW) 

aaa =. 7 
fyi (0 
Both derivatives exist for a smooth system. Using the implicit function theo- 
rem it follows that, provided that the Jacobian matrix J(, ju) is non-singular, 
locally (1) is equivalent to writing x as a function of y, ie, O = f(a(p), 1). 
Then it follows from differentiating (1) with respect to uw that 


dx 


He. w) FE + ful) =0. (8) 


As J(, 1) is non-singular, we can solve for a A point (2, 2) is called regular 
if det(J(ax, w)) #0. 


Definition 1 (Bifurcation). The appearance of a topologically nonequiv- 
alent phase portrait under a variation of parameters is called a bifurca- 
tion [1,6,8, 13, 25]. 


Definition 2. Transcritical bifurcation is a particular kind of local bifurca- 
tion, meaning that it is characterized by an equilibrium having an eigenvalue 
whose real part passes through zero. 

A transcritical bifurcation is one in which a fixed point exists for all val- 
ues of a parameter and is never destroyed. However, such a fixed point 
interchanges its stability region with instability region as the parameter is 
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varied. In other words, both before and after the bifurcation, there is one 
unstable and one stable fixed point [1,6,8, 13]. 


Theorem 1. (Sotomayor Theorem) Suppose that fu (to) = 0 and that n x 
n matric A = Df (xo, Uo) has a simple eigenvalue AX = 0 with eigenvector 
v and that AT has an eigenvector w corresponding to the eigenvalue \ = 
0. Furthermore, suppose that A has k eigenvalues with a negative real part 
and (n —k —1) eigenvalues with a positive real part and that the following 
conditions are satisfied: 


w" fr(xo, Ho) #0, w" [D? f (xo, po(v,v))] #0 (9) 


then there is a smooth curve of equilibrium points for « = f(a,j) in R” x 
Rpassing through (xo, Uo) and tangent to the hyperplane R” x jup. Depending 
on the signs of the expressions in (6). In this case the system experiences a 
saddle node bifurcation. If the conditions (9) are changed to : 


w" f(xo, Ho) = 0, 
w" [Dfu(xo, Ho)] # 0, (10) 
w? [D? f (x0, Ho) (v,v)] # 0, 


then the system (5) experiences a Pitchfork bifurcation. And if the condition 
(9) changed to: 


w" fu(ao, Ho) = 0, w" [D f(x, Ho)v] #0, (11) 
w" [D? f (x0, Ho)(v,v)] = 0, w? [D? f (xo, H0)] 4 0, 


then the system (5) experiences a Transcritical bifurcation. 


Proof. For the proof, one can see [8]. 


3 Bifurcation of the system (12) 
Consider the following system of differential equations: 
(12) 


y=ry(l— %) — ay — pyz 
2 CZ 
2 = id — aye — bz 


where y and z are as before. The virus population is assumed to grow lo- 
gistically: r is the viral replication rate at low viral loads, and we assume 
that this rate is decreased linearly with increased viral load to reach zero at 
a viral load k. Immune cells are assumed to be inhibited by the virus at a 
rate qyz and die at a rate b. 

Clearly Eo = (0,0) is a trivial equilibrium of the system. There exist an 
equilibrium E; = (y,0) = (£(r — a),0) provided r > a > 0. 
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The equilibrium F; is called the virus dominant equilibrium (VDE). More- 
over, we can find another equilibrium E* = (y*, z*), where y* > 0 and z* > 0, 
satisfying the following equations: 


eae 


Wet (13) 


cy 
1+dy* 
E* > 0 means that while the virus population is growing, immune cells 


start to increase; therefore, our main attention will be on equilibrium E*. It 
follows from the first equation of (13) 


f= r(k — y") ~ ak +0 (14) 
pk 


By z* > 0, one can find ¥ such that 
yo <¥ (15) 
In order to find y* for E*, we should solve the quadratic equation 
h(y) = qdy* + (—c+q+bdjy+b=0, y* <y (16) 
to obtain a double root for (16), one should have 


A=0= (c—q-— bd)? = 4bqd > c—q— bd = £2,/bqd (17) 


The minus sign for the root is not applicable, so 
c—q—bd = 2v/bqd (18) 


or equivalently ¢ = (,/¢ + Vbd)?. 
Conditions (17) and (18) on polynomial (16) lead to 


c—q-— bd, 2/bqd .5 b v5 
— = — 1 
av) = (w—- Fo = - EP GP 9) 
Consequently, 
b 
* = 4/— 20 
y ad (20) 
and 
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Because y* < §, we can define a threshold (see the following definition) 
as follow: 


rk—ry* > ak > r(k—y") > ak > —— ; (22) 
—Y 
aks if * k 
>rT = k-y pad = (23) 
oe) ify* >k 


Definition 3. In mathematical or statistical modeling, a threshold model 
is any model where a threshold value, or set of threshold values, is used to 
distinguish ranges of values where the behavior predicted by the model varies 
in some important way. 

With the above statements, one can have the following lemma: 


Lemma 1.Suppose that (18) is satisfied. 

(a) If r <a, then the trivial equilibrium Eo = (0,0) is the only equilibrium. 
(b) Ifa<r<nrni(tie a<rand y* > %), then there are two equilibria Eo 
and E, = (9,0), where y = *(r —a) 

(c) Ifr,<rtiea<r and y* <¥%), then there are three equilibria, Eo, Ey 
r(k—y*)—ak 


and additional equilibrium E* = (y*,2*) with 2* = aR 


We call E* the immune control equilibrium(ICE). 


Here we would like to determine the type of the equilibria (Eo, E; and E*) 
for the system (12). 


3.1 Global dynamics of (12) 


Let (y*, z*) be an equilibrium of (12). The associated characteristic equation 
of (12) is given by 


go(A) = x + cA +Co = 0 (24) 
where 
= ( zr : _ y - b) (25) 
= r a z 
Cl ie y Pp 1 ly* qy 


and 
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2r cy* cz* 


a ee Pe) aie qy” — b) Ser ay 


(26) 


At Eo = (0,0), two roots of the characteristic equation are A; = —b < 0 and 
Ag = —(a—r). Therefore, Eo is stable if r < a. Otherwise, if r > a, then Eo 
is a saddle point. 

At Fy, y* =y, 2* =0, a direct calculation implies that FE is stable. 

At E* we have 


>0 (27) 


and 


SK ( c = * ok 
co = py" 2"( 9+ Tyra? py*z*gi(y)- (28) 


If gi(y) = 0, then g = Vee Substituting g in g(y) where 


d/q 
c—q-—bd\., 
= (y — ———_ 2 
gy) = (y Sad ) (29) 
we have 
P ve—J/q c—q-—bd 
a) = 
qd q 
Vat vod — /g- Vbd,, 
Vd 
(30) 
y = y”, therefore co = 0 
In this case since Ay + Az = —c, < 0, AyAg = 0 which gives us 
1 =0 (31) 
and 
r b 
2 SHS -qy = =i qd (32) 


The system (12) under condition (18) for the equilibrium E* has one nega- 
tive eigenvalue and one zero eigenvalue. Next we check the conditions for the 
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Transcritical bifurcation. For this purpose, we use Sotomayor theorem (see 
Theorem 2.3). In the following, we calculate the Jacobian matrix, second 
derivative of the Jacobian matrix and also eigenvector v corresponding to 
eigenvalue A; = 0 for A and w,the eigenvector of A; = 0, corresponding to 
aX, 

The Jacobian matrix of (12) is 


(tay 9% Tray 8 
where condition (18) implies that 
Cz cy 
(1 + dy)? ‘al 1+ dy Y a) 
Therefore, 
2r 
ee " ey a— pz | (35) 


The Jacobian matrix A at E* will be 


an, Hb be 
Ape = kV ad PY aa (36) 
0 0 


By a direct calculation, the eigenvectors v and w are 


kp 

v= wm) = (21 (37) 
Ww = (wi, W2) = (0,1) (38) 

v1 * 2 * 2 * 2 * 
D? f(B*)(v,v) 2 Ae din, + 2 fe yy + 2 pale Jag, + 2 De dygvy B 

VV) >= 2 * 2 * 2 * 2 * = 

7 Bip ym +2 Be dyjve + 2 Be dyg, + 2 Be ) poVv9 © 
(39) 


If c £0 implies that D? f(E*)(v,v) 4 0. Also, one should have 


b 


ite) =( faa-) (40) 
0 


From (36) and (38), one can obtain 
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T ey 
w f,-(E*) =0 (41) 


The above calculations and results lead to the conclusion that conditions (8) 
are valid. Therefore, by the Sotomayor theorem, the system (10) undergoes 
transcritical bifurcation. 


4 Example (numerical simulation) 


The parameters data are choosen such that the Figures 1-5 are in consistent 
with [2,4,5,9,12,14,15]. Since we are dealing with AIDS, the following Figures 
show the regions of weak immune response. 

We try to find a region for parameter r with respect to parameter k. From 


(194 > ae but y* = ap therefore 


8 (42) 
k—-,/2 
qd 
sO 
b 
k — 43 
> ad (43) 


Thus, the parameter region is obtained in Figure 1: 


Figure 1: Parameter region r with respect to k by considering a = 3, b = 2, q = 9 and 
d=2 


We use numerical techniques to determine the system (12) with condition 
(18). 
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In Figure 2, first we obtain the parameter r with respect to y and z. Next 
by considering the values a = p = 3, k = 4/3, ¢q = 9, b = 2 and r = 10, the 
stability regions of the orbits are investigated. Therefore, system (12) is in 
a steady state; this means that however the immune response of the body is 


so weak that is still can defend against the disease. 


Figure 3 shows that after 100 days, immune cells could not control the growth 


of virus cells and so the system (12) is unstable. 
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Figure 3: a=p=3,k =10, ¢q=9, b=2 and r =5 with initial condition (4, 2) 
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Results Figure 1 shows the region which one can choose r with respect to k. 
In Figure 2, k is small, so one can see the stable regions. Figures 3, 4 and 5 
show an unstable region corresponding to an increase in the virus (k > 10). 
Compares: In this study by assuming condition (4), for system (12), we paid 
attention to AIDS. According to our knowledge, this condition that causes 
a complex dynamic is not considered in any related previous works [2,9, 12]. 
The new condition was lead to an one zero eigenvalue and as a consequence, 
by applying Sotomayor theorem, to transcritical bifurcation. Therefore, we 
determined the stable and unstable regions (by different given values, small 
and large for k) by using transcritical bifurcation. 


5 Conclusion 


In this paper, we analyzed system (12), with condition (4), at the equilibrium 
corresponding to only one zero eigenvalue (co-dimension one bifurcation). 
In order to determine transcritical bifurcation, we applied condition (10) in 
Sotomayer theorem (see theorem 2.3). One can notice that as we mentioned 
in the results of our investigation, the difference between this study and 
others [2,9, 12]. From the biological point of view, the stable and unstable 
regions correspond to the viral population load. Moreover Figures 3, 4 and 
5 showed that the virus population of AIDS increases for the value of k > 10. 
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